## 

combine_imputation_models = function(reg.list){
  # reg.list is a list of regression models, assumed to be from runs of imputation 
  coef.mat = se.mat = NULL
  
  for(reg in reg.list){
    
    cm = summary(reg)$coefficients
    coef.mat = cbind(coef.mat, cm[,1])
    se.mat = cbind(se.mat, cm[,2])
    
  }
  
  mean.estimates = apply(coef.mat, 1, mean)
  if(ncol(coef.mat) == 1){
    bw.var = 0
  }else{
    bw.var = apply(coef.mat, 1, var)*(1 + 1/ncol(coef.mat))
  }
  avg.wi.var = apply(se.mat^2, 1, mean)
  se.estimates = sqrt(avg.wi.var + bw.var)
  df = data.frame("term" = rownames(cm), "estimate" = mean.estimates, "std.error" = se.estimates)
  rownames(df) = NULL
  
  df
  
}